Assessment

#How did covid-19 pandemic and lockdown affect the prescription trends of the most common antibiotics in Scotland.

I chose this question because I wanted to look into the effect of COVID-19 on the prescription of drugs across Scotland. I chose to analyse the antibiotic trends because I hypothesised that covid 19 would cause a decrease in the prescription of antibiotics due to the focus on COVID-19 treatment and the lifestyle during the main 2020 lockdown. This lifestyle included mask-wearing, social distancing and the use of hand sanitiser. These changes would have caused a reduction in the spread of both COVID-19 and bacteria, therefore reducing the prescription of antibiotics. I also hypothesise that by 2024, the number of prescriptions will have returned to a similar level in 2019. This is because doctors may want to keep antibiotic prescriptions low to prevent antibiotic resistance, but this effect is cancelled out by the growing population.

One scientific journal reinforced this theory by analysing data from Poland. They found that during COVID-19 there was a decrease in the number of medical consultations (both online and in-person) and a decrease in total prescription of antibiotics. Similarly, in Portugal they found a decrease in antibiotic prescriptions during COVID-19 (Silva et al., 2021). In England, it was found there was a 12.4% decrease in the prescription of antibiotics for respiratory tract infections in COVID (Andrews et al., 2021).

I used the post-pandemic data to be from 2024 and not 2021, because I wanted to make sure that there were minimal COVID illnesses or restrictions that could affect the results.

To answer this question, I first needed to load the appropriate libraries. Then, I loaded the datasets and created a function to concisely filter and edit the datasets for each year. I learned how to do this from the website https://www.dataquest.io/blog/write-functions-in-r/. In this function, I filtered for antibiotics by only looking at prescriptions with the BNFItemCode starting with 0501 (Open Prescribing, 2019). Then I filtered out NAs, cleaned up the names and calculated the number of prescriptions and arranged the data in descending order.

#Load the data for each year I will be looking at
may2020 <- read_csv(here("data","pitc202005.csv"))
may2024 <- read_csv(here("data","pitc202405.csv"))
may2019 <- read_csv(here("data","pitc201905.csv"))

#To make the relevant column titles consistent across datasets
may2019 <- may2019 %>% 
  rename("HBT" = HBT2014)

#A function was created to filter and edit the datasets. This prepares them for joining
process_antibiotics <- function(data) {
  data %>%
    filter(str_detect(BNFItemCode, "^0501"), !is.na(BNFItemDescription)) %>% 
    mutate(BNFItemDescription = str_replace(BNFItemDescription, "_", " "),
           Antibiotic = str_extract(BNFItemDescription, "^[^ ]+")) %>% #Extracts the first word, which is the antibiotic name
    group_by(Antibiotic) %>%
    summarise(quantity_sum = sum(PaidQuantity), .groups = "drop") %>%
    arrange(desc(quantity_sum)) %>% 
    mutate(quantity_sum = round(quantity_sum))
}

#Apply the function
may2019_processed <- process_antibiotics(may2019)
may2020_processed <- process_antibiotics(may2020)
may2024_processed <- process_antibiotics(may2024)

I prepared the data for joining by renaming Fluclox to Flucloxacillin. I then joined the data, and tidied it. I sliced the data so I was focusing on the top 10 most prescribed antibiotics in 2019 in comparison to 2020 and 2024.

#Changing Fluclox to Flucloxacillin
may2020_processed <- may2020_processed %>% 
  mutate(Antibiotic = ifelse(Antibiotic == "FLUCLOX", "FLUCLOXACILLIN", Antibiotic))

may2019_processed <- may2019_processed %>% 
  mutate(Antibiotic = ifelse(Antibiotic == "FLUCLOX", "FLUCLOXACILLIN", Antibiotic))

#Joining the datasets
antibiotics1920 <- full_join(may2019_processed, may2020_processed, by = "Antibiotic")

antibiotics192024 <- full_join(antibiotics1920, may2024_processed, by = "Antibiotic")

#Each year was renamed and slicing to compare the top 10 most common antibiotics in 2019 to their prescription numbers in 2020 and 2024
antibiotics192024 <- antibiotics192024 %>% 
  rename(
   "May 2019" = quantity_sum.x,
   "May 2020" = quantity_sum.y,
   "May 2024" = quantity_sum) %>% 
  slice(1:10)

#The antibiotic names were changed to only have the first letter capitalised. I learnt this from the R help page
antibiotics192024 <- antibiotics192024 %>% 
  mutate(Antibiotic = str_to_sentence(Antibiotic))

I created a stacked bar graph with ggplot and plotly to visually compare the prescriptions.

options(scipen=999)

#The antibiotics data was pivoted to add an extra column for the year. This is so there can be separate bar graphs for each year
pivoted_antibiotics <- antibiotics192024 %>% 
  pivot_longer(
    cols = c("May 2019", "May 2020", "May 2024"), 
    names_to = "Year", 
    values_to = "Prescription_Number")

#Creating the graph from the pivoted data
antibiotics_graph <- pivoted_antibiotics %>% 
  ggplot(aes(x = Prescription_Number, y = Year, fill = Antibiotic)) +
  geom_bar(stat = "identity") +
  labs(title = "Prescriptions Sold for Each Antibiotic by Year") +
  scale_fill_brewer(palette = "Paired") +
  scale_x_continuous(labels = comma) +
  theme_minimal()

#ggplotly was used so when a mouse hovers over the bar graph, the prescription name and number can be seen
antibiotics_graph <- antibiotics_graph %>% 
  ggplotly(tooltip = c("Prescription_Number", "Antibiotic"))
antibiotics_graph

This data visualisation shows that the overall number of prescriptions sold decreased during lockdown and was slightly lower before lockdown than after. It also shows that some antibiotics had a larger reduction in prescriptions than others. For example, Trimethoprim, Lymecycline and Flucloxacillin seem to stay somewhat consistent in 2019, 2020 and 2024. However, Phenoxymethylpenicillin and Amoxicillin seem to fluctuate more. This could be because different antibiotics treat different infections. Trimethoprim treats UTIs (Urinary Tract Infections) (NHS, 2019), so it is likely it stayed consistent through COVID-19 because the bacteria that caused the UTI came from the individual, and so lockdown didn’t affect its spread. While, Phenoxymethylpenicillin is a type of Penicillin that treats upper respiratory infections, which are contagious (NHS, 2022). This explains why the lockdown would’ve affected the need to prescribe some types of antibiotics and not others.

I then created a gt() table to show the number of prescriptions of each antibiotic in each year, and calculated the total prescriptions per year.

#To find the total antibiotics across each year and bind it to the dataframe
all_antibiotics_table <- antibiotics192024 %>%
  summarise(across(starts_with("May"), sum, na.rm = TRUE)) %>%  
  mutate(Antibiotic = "Total") %>%
  bind_rows(antibiotics192024, .)

#To create a gt() table
all_antibiotics_table <- all_antibiotics_table %>%
  gt() %>% 
  cols_label(Antibiotic = "Antibiotic") %>% 
  cols_align(align = "center",
             columns = everything()) %>% 
  tab_style(
    style = cell_text(weight = "bold"),
    locations = cells_column_labels(columns = everything())
  ) %>% 
  tab_header(
    title = "Affect of the Covid19 pandemic on the prescription of 2019's top 10 most common antibiotics",
    subtitle = "Looking at how the number of prescriptions for the top 10 most common antibiotics in 2019 changed during and after the 2020 Covid19 lockdown") %>% 
  tab_spanner(
    label = "Number of Prescriptions",
    columns = c("May 2019", "May 2020", "May 2024")) %>% 
  fmt_number(
    columns = c("May 2019", "May 2020", "May 2024"),
    decimals = 0,
    use_seps = TRUE  # To add commas to make the number more easily readable
  )
all_antibiotics_table

Affect of the Covid19 pandemic on the prescription of 2019's top 10 most common antibiotics
Looking at how the number of prescriptions for the top 10 most common antibiotics in 2019 changed during and after the 2020 Covid19 lockdown
Antibiotic	Number of Prescriptions
	May 2019	May 2020	May 2024
Amoxicillin	2,875,868	1,778,491	3,751,932
Phenoxymethylpenicillin	1,892,108	968,754	2,616,743
Flucloxacillin	1,612,959	1,475,567	1,733,244
Trimethoprim	782,422	745,904	696,670
Erythromycin	661,545	443,204	582,068
Lymecycline	623,567	493,844	538,300
Oxytetracycline	572,089	428,371	277,441
Doxycycline	447,389	387,163	616,170
Co-amoxiclav	382,770	349,958	466,694
Clarithromycin	377,144	230,282	493,464
Total	10,227,861	7,301,538	11,772,726

This gt() table illustrates the number of prescriptions for each antibiotic and the total prescriptions in 2019, 2020 and 2024. This table helps indicate how much lower the prescriptions were in 2020 and helps directly compare the prescriptions of each antibiotic between years. This table also reinforces the conclusions made from the stacked bar graph, showing that in 2020 there was a clear decrease in prescriptions and 2024 had more prescriptions than 2019 and 2020.

To prepare for creating a geospatial map the healthboard data was loaded, prepared for joining, cleaned, organised and joined to the antibiotics dataframe.

#Load data
NHS_healthboards <- st_read(here("data", "SG_NHS_HealthBoards_2019.shp"), quiet = TRUE)

#Prep for joining
NHS_healthboards <- NHS_healthboards %>% 
  rename("HBT" = HBCode)

#Summarise antibiotic data by HBT code to give a total number of prescriptions per area
antibiotics_by_HBT <- function(data_frame) {
    data_frame %>% 
      group_by(HBT) %>%
      summarise(quantity_sum = sum(PaidQuantity), .groups = "drop") %>%
      arrange(desc(quantity_sum))
}

#Apply the function
HBTmay2019_processed <- antibiotics_by_HBT(may2019)
HBTmay2020_processed <- antibiotics_by_HBT(may2020)
HBTmay2024_processed <- antibiotics_by_HBT(may2024)

#Join antibiotic data
HBTantibiotics1920 <- full_join(HBTmay2019_processed, HBTmay2020_processed, by = "HBT")
HBTantibiotics192024 <- full_join(HBTantibiotics1920, HBTmay2024_processed, by = "HBT")

#Rename and remove NAs
HBTantibiotics192024 <- HBTantibiotics192024 %>% 
  rename(
   "May 2019" = quantity_sum.x,
   "May 2020" = quantity_sum.y,
   "May 2024" = quantity_sum) %>% 
  filter(!if_any(everything(), is.na))

#Pivot to make faceting by year easier
HBTpivoted_antibiotics <- HBTantibiotics192024 %>% 
  pivot_longer(
    cols = c("May 2019", "May 2020", "May 2024"),
    names_to = "Year", 
    values_to = "Prescription_Number")

#Join with healthboard data
Healthboards_and_antibiotics <- NHS_healthboards %>% 
  full_join(HBTpivoted_antibiotics, by = "HBT")  

I then created a faceted geospatial map to show which areas of Scotland were prescribed more antibiotics than other areas, and how this changed during and after Covid19

#To create the geom plot graph
Healthboards_and_antibiotics_plot <- Healthboards_and_antibiotics %>% 
  ggplot(aes(fill = Prescription_Number)) +
  geom_sf(color = "white", linewidth = 0.1) +
  scale_fill_viridis_c(name = "Total No. of Prescriptions",
                       labels = scales::comma) +
  labs(title = "Total Number of Antibiotic Prescriptions in Areas of Scotland from 2019 to 2024") +
  facet_wrap(~ Year) +
  #This code was to create a north arrow and scale on the plot
  annotation_scale(    
    location = "tl"
  ) +
  annotation_north_arrow(
    location = "tl",    
    pad_y = unit(0.5, "in"),    
    style = north_arrow_nautical(
      fill = c("grey40", "white"),      
      line_col = "grey20"
    ) 
  )
Healthboards_and_antibiotics_plot

This geospatial map was used to compare the prescription trends across different areas of Scotland before, during and after the COVID-19 pandemic and resulting lockdowns. This map illustrates that there was an overall decrease in the prescription of antibiotics across Scotland. Notably, Greater Glasgow and Clyde and Lothian have had a large increase in antibiotic prescriptions in 2024. This could be because the big cities of Edinburgh and Glasgow are found in these areas, so even if there is a much larger increase in prescriptions in these areas, it is likely they have proportionally a similar increase in prescriptions compared to other less populated regions like the highlands.

The data shows that my hypothesis that the prescriptions of antibiotics during COVID-19 decreased is true, but my hypothesis that prescription levels in 2024 are similar to those in 2019 is false. This hypothesis was false because the prescription levels in 2024 seem to be higher than in 2019. This increase may be due to COVID-related secondary bacterial infections. Shafran et al. (2021) found that COVID-19 patients had more secondary bacterial infections than those with influenza.

To further assess the data, the analysis of antibiotic prescriptions from 2020 to 2024 could be analysed to see if there was a steady or harsh increase in antibiotic prescriptions after the pandemic. Furthermore, the data could be analysed while taking population into account, to make it clear in geospatial maps if there are some areas that are disproportionately being prescribed more or fewer antibiotics during the pandemic. This data could be used to assess if antibiotics are currently being overprescribed, and could therefore prevent further antibiotic resistance. In the long term, this data could be continuously compared to incoming data for 2025, 2026, etc. to see if there are any further long-term effects of the COVID-19 pandemic on the prescription of antibiotics.

To conclude, my hypothesis that COVID-19 caused a decrease in antibiotic prescriptions is true, however I was incorrect in hypothesising that in 2024 there would be a somewhat equal level or antibiotic prescriptions. For short-term research, the data could be analysed to look at region-specific data in the context of the population. For long-term research, the current data could be compared to future data to see the long-term effects of COVID-19 on prescriptions of antibiotics.